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ABSTRACT OF THE THESIS 


Simulated Groundwater Tracer Study of the Alamitos Barrier 
Project, Los Angeles County, California 

by 

Joseph Christopher Pope 

Master of Science in Civil Engineering 
University of California, Los Angeles. 2006 
Professor William W-G. Yeh, Chair 

A simulated tracer study is performed on the Alamitos Groundwater Barrier 
Project (ABP), Los Angeles. California. The ABP consists of 43 injection wells that 
create a freshwater barrier to prevent intrusion of seawater into the coastal aquifer of 

I 

southern Los Angeles County. The injection wells currently use a blend of recycled 
and imported water. Regional water quality regulations dictate the quantity of 
recycled water that can be injected into the aquifer. The regulations also require that 
the travel-time between the injection wells and the local production wells be greater 
than one year. The purpose of this study is to analyze the travel time of injected 
water into the ABP. A previously calibrated three-dimensional, finite element, 
coupled groundwater flow and transport model is used to simulate the movement of a 
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conservative tracer in the Alamitos Barrier Project. The results of the simulations 
show that tracer travel times between the injection wells and the production wells 
typically exceed 60 years. The model results exhibit a high amount of numerical 
error when large time steps are used. Time-step sensitivity analysis indicates that a 
time step of one day or less will minimize model error. Further analysis of local head 
boundary conditions is recommended in order to provide a wider range of tracer 
travel times given variations in boundary hydraulic head levels. 


IX 


I INTRODUCTION 


Since the late 19"’ Century, groundwater has been an important source of 
water for the Southern California region. From the 1940’s through the 1950’s a rise 
in chloride levels was noticed in many coastal pumping wells. The rise in chloride 
levels was due to the intrusion of seawater into the coastal aquifers as a result of the 
significant depletion of fresh groundwater (CalUson et al. 1991). In the late 1950’s 
and early 1960’s the Los Angeles County Department of Public Works (LACDPW) 
constructed and began operating three seawater barriers along the County’s coastline: 
the West Coast Basin; Dominguez Gap; and Alamitos Barrier Project (ABP). The 
three barrier projects serve to protect inland groundwater reservoirs from 
contamination due to seawater intrusion and also act as a source of recharge water for 
the area aquifers. The ABP, the subject of this study, has been in operation since 
1964 and consists of 43 injection wells and four extraction wells. The injection wells 
create a freshw ater pressure ridge that prevents further intrusion of seawater, while 
the extraction wells create a trough between the coastline and the injection wells, 
further amplifying the effect of the freshwater pressure ridge (Water Replenishment 
District of Southern California 2006). The ABP is located between the southern 
coastal boundar>' of Los Angeles County and the northern coastal boundary of Orange 
County. An overview of the ABP location is shown in Figure 1. 
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Figure 1. Location of Alamitos Barrier Project 

(Water Replenishment District of Southern California 2006) 


The Water Replenishment District (WRD) of Southern California and the 
Orange County Water District (OCWD) supply the injection water used for the ABP. 
The WRD provides injection water for wells on the Los Angeles County side of the 
project, while OWCD supplies injection water to the Orange County wells. The 
water used for injection is purchased from the Metropolitan Water District of 
Southern California (MWD) and is typically a blend of Colorado River and State 
Water Project water. During Water Year (WY) 2004-2005,4,555 acre-feet (AF) of 
imported water was injected into the ABP. Up until WY 2004-2005 only imported 
water was injected into the ABP. However, WRD began augmenting this supply with 
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recycled water in October 2005 (WY 2005-2006). All recycled water injected into 
the ABP is treated to meet or exceed drinking water standards prior to injection 
(Water Replenishment District of Southern California 2006). In order to satisfy the 
permit requirements of Regional Water Quality Control Board (RWQCB) and the 
Slate Department of Health Services (DHS), WRD is required to provide an analysis 
of the temporal and spatial distribution of injection water in the ABP. The standards 
set forth by the RWQCB state that no more than 50% of the water injected into the 
aquifer can be from recycled water and the travel time of the injected water between 
the barrier and any production well must be at least one year. 

In addition to the 43 injection and four extraction wells that make up the 
barrier project, there are four production wells in the study area that provide water for 
primarily industrial uses. The main objective of this study is to analyze the travel 
time of injected water from the barrier injection wells to the four production wells. 
This study utilizes the results from the dissertation by Bray (2006). Bray simulated 
the physical system of the ABP using FEM WATER version 3 (Lin et al. 1997), an 
open source 3-dimensional, finite element, coupled groundwater flow and transport 
model. The model was calibrated using borehole data and 11 years of historical head 
and chloride concentration data taken from a network of 180 observation wells, along 
with the ABP historic injection and extraction rates. The study herein uses the 
calibrated flow and transport model to measure the travel time of an injected 
conserv'ative tracer, which is postulated to be representative of the movement of the 
reclaimed water in the ABP. 
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2 LITERATURE REVIEW 


There are a multitude of numerical techniques available today that give water 
resource managers the ability to characterize the flow of groundwater in geologically 
complex, multi-layered aquifers. The reasons for characterizing the flow in a given 
groundwater basin are numerous, including protection of groundwater supplies from 
environmental degradation and over-pumping. A study of the groundwater flow by 
Cunningham et al. (1994) at Wright-Patterson Air Force Base, Ohio was conducted 
using the particle tracking module (MODPATH) of the United States Geological 
Survey’s modular three-dimensional finite-difference ground-water flow model 
developed by McDonald and Harbaugh (1988) (hereafter referred to as 
MODFLOW), The results of the study by Cunningham et al. (1994) illustrated the 
groundwater path lines as they crossed the Base, giving Air Force planners a tool that 
could assist them in developing future remediation plans. Snyder et al. (1996) also 
utilized the MODPATH module of MODFLOW in conducting a groundwater 
vulnerability study for Clark County, Washington. The results of their study was 
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useful in identifying parts of the groundwater flow system that might be affected by 
effluent from on-site waste-disposal systems, in addition to estimating the age of 
groundwater for any part of the system. 

Tompson et al. (1999) analyzed groundwater residence times in a large 
groundwater recharge operation in Orange County, California. Using a 2- 
dimensional finite difference model, their study showed that groundwater residence 
times were significantly shorter than residence times inferred using naturally 


4 


occurring isotopic tracers. The large difference in residence times between the two 
schemes was due to the numerical model’s incomplete representation of the geologic 
heterogeneity of the groundwater system. The groundwater system modeled in the 
Orange County study is similar to the ABP with regards to both geomorphology and 
scale. 

Daniel et al. (2000) used a multi-scale statistical framework to estimate 
groundwater travel times and to derive travel time and probability densities. They 
showed that travel time uncertainties depend primarily on uncertainties in hydraulic 
conductivity and that the probability densities were best approximated using a Monte 
Carlo technique. Wen and Kung (1996) used the constant displacement scheme in the 
random walk method as a means of eliminating numerical dispersion and increasing 
computational efficiency for particle tracking in a heterogeneous aquifer. Lodiciga 
(2004) modeled the total travel distance and hydraulic conductivity in steady-state 
groundwater flow systems using asymmetric gamma distribution functions. 

In the current study the problem of groundwater travel time is solved using the 
coupled groundwater flow and transport equations for a single constituent. One of the 
common problems in transport modeling is the propagation of errors due to numerical 
dispersion. Numerical dispersion can be reduced by using a variety of techniques. 
Satisfaction of the Courant condition helps to minimize numerical dispersion in 
transport problems. The Courant condition for one-dimensional steady state flow is 
shown in Equation 1. 
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Equation 1 - Courant Condition, 

where At is the time step, Ax is the spatial discretization in one-dimension, and v is 
the average velocity in the x-direction. The Courant condition ensures that a particle 
of water does not travel further than Ax in one time step. Satisfying the Courant 
condition can be challenging in an unsteady, complex 3-dimensional finite element 
model, where the discretization of the physical model domain can vary from very fine 
around points of interest, such as wells, to relatively coarse discretization near the 
boundaries. Woods et al. (2003) analyzed the numerical error associated with the 
finite element transport code, SUTRA (Voss 1984). They developed a method for 
subtracting numerical dispersion from the SUTRA solution, but found that for 
complex nonlinear transport problems SUTRA produced an unacceptable amount of 
numerical error even with a very small At and a highly refined grid. For the sake of 
computational efficiency others have developed solution algorithms that take into 
account the unsteady nature of the flow solution and the spatial variability of the 
discretized model space (Bensabat 2000; Binning and Celia 2002; Huang 1992; 
Kaasschieter 1995; Zheng 1994). 

Clarke et al. (2005) developed a probabilistic, one-dimensional travel 
distance/time equation, known as the “fractional advection-dispersion equation” 
(ADE). The fractional ADE is based on the concentration that one wishes to know 
the travel distance or time for, and the value of a parameter chosen based on the 
degree of heterogeneit)' of the aquifer. Clarke et al. found that the fractional ADE. 
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provided better estimates of travel distance or time than for estimates based on 
traditional numerical solutions of Darcy’s law. 

3 METHODOLOGY 

3.1 Conceptual Model Description 

The study area encompassing the ABP is approximately 11 square miles. The 
Seal Beach and Los Alamitos Faults represent the southwest and northeast boundaries 
of the model, respectively. The model boundaries are extended several thousMid feet 
away from the barrier in order to reduce boundary effects. A plan view of the model 
area with the overlaid finite-element mesh is shown in Figure 2. The production 
wells vary between one and two miles from the injection wells. 

There are 11 geologic layers represented by 23 numerical layers, starting with 
an upper unconfined aquifer followed by a series of aquitards and aquifers. Both the 
geologic layers and the numerical representation of the model layers are shown in 
Figure 3. 
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Figure 2. Alamitos Barrier Project Map 
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Figure 3, A BP Geologic Layers 
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The original conceptual model developed by Bray (2006) included only the 
layers directly affected by the injection wells and the upper unconfined aquifer. For 
the purposes of this study the Main Zone, or lower aquifer, was added to the 
conceptual model. Outside of the mode! domain, the Main Zone aquifer is a 
productive groundwater source for the local municipalities, therefore it is important to 
analyze the impact of injection on this part of the system. Due to a limited amount of 
data pertaining to the geometry and geophysical characteristics of the Main Zone 
aquifer and Main-I aquitard. these layers are represented in the model as a loose 
approximation of the I aquifer and I-A aquitard. The ABP model is represented in 
FEMWATER as a finite element mesh containing 120,000 nodes that define 227,723 
finite elements. 

3.2 Groundwater Flow and Transport Simulation 

The FEMWATER model developed by Lin et al. (1997) utilizes a density- 
dependent flow and transport solution and is represented as follows: 



Equation 2 - Density-dependent flow model 



+ F,.«VC-V* 
di 



Equation 3 - Advection-dispersion model 
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Equation 4 - Storage coefficient 


M 

Equation 5 - Hydraulic conductivity tensor 


OD = a, V S + (oT/ - a, )-j^ + a„0TS 


Equation 6 - Moisture Content / Dispersion Relationship 


f/, + 

Equation 7 • Darcy flux equation 


— = a^ +i/,C 
Pa 

Equation 8 - Relationship between concentration and density 

Equation 2 is the conserv'ation of mass equation where h [I,] is related to the 
fluid density p [ML‘^], the fresh water density pn the density of the injected 

fluid p* [ML‘^], a source or sink q [T’j, and time / [Tj. The primary parameters in 

Equation 2, the specific storage F [L‘‘] and the hydraulic conductivity tensor K [LT~ 
'], are calculated using Equations 4 and 5. The a coefficient in Equation 4 is the 
coefficient of consolidation of the soil [LT^M''], g is the gravitational acceleration 
constant, $ is the moisture content, n is porosity, ^ is the compressibility of the fluid 
[LT^M '], and the final term describes the differential change in saturation {S) with 

head. In Equation 5,p [ML''T~'] is the dynamic viscosity of water and k [L^] is the 
intrinsic permeability tensor. 
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Equation 3 shows the advection dispersion model where C [ML*^] is the 
constituent concentration, m is the external source or sink rate [ML'^T'], Vp [LP’] 
is the Darcy velocity (computed in Equation 7), and D [L^P'l is the dispersion 

I 

tensor. Equation 6 is used to calculate the third term in Equation 3, which is the 
gradient of the surface flux, where ai [L] is the longitudinal dispersivity, aj [L] is the 
transverse dispersivity, S is the Dirac delta fiinction, V [LP'] is the fluid velocity and 
is proportional to Pp, a„ [L^P‘] is the molecular diffusion coefficient, and r is the 
tortuosity. The molecular diffusion coefficient is typically assumed to be negligible. 

Equation 8 is the first-order relationship that links the variation in 
concentration to fluid density. For the purposes of this simulation, the injection water 
is treated as a conservative tracer with the same density as the fresh groundwater, 
therefore a/ is set to unity and is set to zero, which effectively decouples the flow 
and transport models in FEMWATER. Additionally, adsorption of the tracer to the 
aquifer medium is assumed to be zero. A list of the spatial and time invariant 
parameters utilized in the simulation is shown in Table 1. 
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Table 1 - List of Spatial and Time Invariant Parameters 


PARAMETER 

VALUE [Unitsl 

Pfl, freshwater density 

1.9383 [slugs/ft^] 

p, dynamic viscosity of water 

2.347 [slugs/ft/d] 

g, acceleration of gravity 

32.2 [ft/s"] 

p, compressibility of water 

2.822x10'”* [ft-d^/slugs] 

a. compressibility of aquifer medium 

3.27x10-'" [ft-d"/slugs] 

n, porosity of aquifer medium 

0.43 [-] 

ol, longitudinal dispersivity 

50.0 [ft] 

Ct, transverse dispersivity 

5.0 [ft] 

T, tortuosity of aquifer medium 

1.0 [-] 


The hydraulic conductivity varies spatially depending on local geologic 
conditions. The hydraulic conductivity field was determined using geologic data 
from 179 boreholes and a natural-neighbor-kriging interpolation scheme. The 
longitudinal and transverse dispersivity values were determined during the model 
calibration by Bray (2006). 

The Dirchlet boundary condition is the prescribed condition for all model 
boundaries. Both the initial and boundary' conditions for head were determined from 
historic observations of head from 1992 (Bray 2006). The initial condition for tracer 
concentration was set to zero everywhere within the model domain. The volumetric 
injection rate at each of the barrier injection wells was set to the maximum injection 
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capacity of each individual well. These values were determined based on the 
maximum historic injection rates. The tracer concentration injected into each well 
was set at 1000 parts per million (ppm). This value was set arbitrarily high in order to 
reduce numerical error. Higher tracer values were experimented with, but the model 
had difficulty converging to a solution at injection values greater than 1000 ppm. 

3.3 Simulation Scenarios 

The initial simulations were run with a one-year time step and ten-year 
duration. After a ten-year simulation the resulting tracer concentrations at the 
production wells were effectively zero, so a longer duration simulation was utilized. 
The primary simulations were then run for 100 years using both one-month and one- 
year time steps. In order to test the model sensitivity to changes in time-step, seven 
runs were performed over a 10-year duration using time-steps that varied from a half¬ 
day to one-year. The concentration breakthrough curve for the time-step sensitivity 
analysis was taken at an observation node set relatively close to the barrier injection 
wells. All simulations were run on a Linux cluster node operating two Intel 2,2 GHz 
Xeon processors with 2 gigabytes of random access memory. It took approximately 
48-hours to run a 100-year simulation using one-month time steps. Post-processing 
was performed using a commercial software package. 


4 RESULTS 


4.1 Time-Step Sensitivity Analysis 

The breakthrough curves for concentration at an observation node located 
approximately 600 feet from the nearest injection well are shown in Figure 4. The 
selection of an observation node relatively close to the barrier was chosen in order to 
minimize the higher computation time required by the smaller time steps. The results 
show that as the time step increases, the observed concentration also increases. The 
results indicate that a one-day time step produces marginally better results than a half¬ 
day time step. Making the assumption that the half-day time-step simulation 
produced an acceptable level of numerical dispersion, an error analysis is performed 
comparing it to the one-day through 360-day time step simulations. The formulation 
of the error analysis is shown below in Equation 9. 



RMSE = 


Equation 9 - RMSE Calculation, 


where n is the number of years, is the "unbiased” concentration for the half-day 

time-step simulation at year /, and C, is the observed concentration at year i for all 
other simulations where time-step is between one day and one year. 
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Graph of Tracer Concentration at Observation Node 80908 over 10 year Simulation Period 
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Figure 4 - Time-Step Size Sensitivit}' Analysis 
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Figure 5 - Root Mean Squared Error Analysis 
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4.2 Simulated Tracer Concentration at Production Wells 

The concentration breakthrough curves for the 100-year simulations using 
time-steps of one month and one year were examined for each production well. 
Figure 6 shows the locations and associated names of each of the four production 
wells in the ABP: SCWC-LAYT, SB-BEV, SB-LEI, and RUIZ-6F1. 


i-. 

Figure 6 - Location of Production Wells 


RUI2-6F1 
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Figures 7 through 10 show the concentration breakthrough curves in parts per 
trillion (PPT) at each of the production wells. For well SCWC-LAYT a 100-year 
simulation shows that for a time step of one-year the tracer concentration begins to 
increase after approximately 60 years. At the end of the simulation the concentration 
is 4.7x10"* PPT, or 47 trillionth of a percent of the injected tracer. For SCWC-LAYT, 


the 100-year simulation using a monthly time-step showed that after 100 years the 
tracer concentration was roughly two orders of magnitude smaller than for the one 
year time-step simulation, 



Figure 7 - Concentration Breakthrough Curve for SCWC-LAYT 

The tracer concentration breakthrough curves for well SB-BEV are shown in 
Figure 8. The yearly time-step simulation for SB-BEV shows that concentration 
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levels begin to increase around the 60-year point and reach about 0.002 PPT, or 0.2 
billionth of a percent of the injected tracer. The monthly time-step simulation for SB- 
BEV shows a nominal rise in concentration levels after 90 years. After 100 years the 
monthly time-step simulation show'ed a concentration of 2.4x10'^ PPT, three orders of 
magnitude lower than for the yearly time-step simulation. 



Figure 8 - Concentration Breakthrough Curve for SB-BEV 

The tracer concentration breakthrough curves for well SB-LEI are shown in 
Figure 9. The yearly time-step showed the appearance of tracer around 55 years and 
a concentration value of approximately 1.3 PPT after 100 years. The tracer 
concentration in SB-LEI after 100 years is 1.3 billionth of a percent of the injected 
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tracer concentration. The monthly time-step simulation shows concentration of tracer 
beginning to rise after 75 years and reaching 0.08 PPT after 100 years. 



Figure 9 — Concentration Breakthrough Curve for SB-LEI 

Figure 10 shows the tracer concentration breakthrough curves for well RUIZ- 

I 

6F1. The yearly time-step simulation shows an increase in concentration in RUIZ- 
6F1 after 85 years with the concentration reaching 0.04 PPT after 100 years. The 
monthly time-step scenario tor RU1Z-6F1 shows concentration levels rising after 
around 65 years and reaching 0.05 PPT after 100 years. The results for this well 
differ from the trends shown in the breakthrough curves of the previous three wells. 
This difference could be due to the close proximity of RU1Z-6F1 to the barrier 
injection wells. For RUIZ-6F1 the monthly time-step produces a slightly more 
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conservative result. The 100-year concentration level is 10% higher for the monthly 
time-step than for the yearly time-step simulation. 



Figure 10 - Concentration Breakthrough Curve for RUIZ-6F1 
4.3 Velocity Field Analysis 

To provide further information on the nature of the groundwater flow in the 
ABP it is insightful to analyze the velocity fields for each aquifer layer. TTie 
groundwater velocities in the entire system varied between 1.3x1 O'* ft/day to 12 
fl/day. The average nodal velocity was 0.6 ft/day and the median was 0.17 ft/day. 
The velocity fields shown in the Figures 11 through 14 are snapshots taken after a 
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one-year simulation. The highest velocities occur near the injection wells for all 
aquifer layers where injection occurs. 

Figure 11 shows the velocity field in the x-y plane for the B-aquifer. Well 
RUIZ-6F1 is the only production well that pumps out of the B-aquifer. From the j 
velocity vectors the flow pattern in the B-aquifer shows that the flow is generally 
outward from the injection wells and then trends towards the north and northwest. 
Well RUIZ-6F1 is approximately 5,000 feet from the nearest barrier injection well. 
From the injection well to the production well, the flow generally turns northward 
before reaching the production well. 
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Figure 12 shows the velocity field in the x-y plane for the A-aquifer. 
Production wells SCWC-LAYT and RU1Z-6F1 are the only two wells that pump 
from the A-aquifer. The flow pattern in the A-aquifer is similar in nature to that seen 
in the B-aquifer where there is outward flow in all directions from the injection wells 
with the flow trending towards the north on the eastern side of the barrier. 



Figure 12 - A-Aquifer Velocity' Field 

Figure 13 shows the flow pattern in the l-aquifer. Three of the four 
production wells pump water from the 1-aquifer. The l-aquifer exhibits a similar flow 
pattern as the A and B aquifers with the exception of a strong gradient from the 
southeast boundary, which causes the flow direction around the production wells to 
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trend strongly towards the northwest. This flow pattern from the southeast boundary 
effectively limits the influence of the injection wells on the production wells. 



Figure 13 - I-Aquifer Velocity Field 

The velocity field in the Main Aquifer is shown in Figure 14. Since there is 
no pumping or injection taking place in this aquifer the velocities are relatively low 
and the flow pattern driven predominantly by boundary conditions and influence from 
the overlying 1-Aquifer. 
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5 DISCUSSION 

As shown in Figures 7 through 10, the travel time of the tracer from the 

I 

injection wells to the production wells varies greatly depending on the size of the time 
step. Using a one-year time step, the tracer travel time varies between 55 and 85 
years. The one-month time step simulation shows that travel time varies between 65 
and 100 years. The large disparity' in travel times between the two simulation 
methods is due primarily to the propagation of numerical error when using the larger, 
annual time-step. The time-step sensitivity analysis supports the conclusion that a 


24 





smaller time-step will produce less numerical error. However, there is a significant 
computational expense associated with a finer temporal discretization. 

Analysis of the velocity fields in each aquifer show that the predominant 
influencing factors on the flow patterns are the boundary conditions and the high 
hydraulic gradients created by the barrier injection wells. The constant boundary 
head condition that was used for all simulations is not necessarily a realistic boundary 
condition for a 100-year simulation. Due to sparse historical data on aquifer head 
levels, determination of the true boundary conditions can be difficult. Further 
analysis utilizing boundary conditions representing drought conditions would provide 
a wider range of potential tracer travel times. The velocity patterns seen in all 
aquifers, the long tracer travel-time, and the extremely low tracer concentrations seen 
at the production wells indicate that the tracer transport is driven predominantly by 
dispersion. Due to the dispersion dominated transport of the tracer, additional 
analyses of the model sensitivity to changes in the dispersivity parameters may be 
useful. Assuming that the injection rates remain relatively constant, dropping the 
head levels at the boundaries would likely increase the advective transport of the 
simulated tracer toward the production wells, which would effectively decrease tracer 
travel times. The constant pumping and injection rates, if varied, may also have an 
effect on the tracer travel times. Reduction in both the injection and pumping rates 
would likely increase tracer travel times. 

Based on the simulations that were performed in this study, the requirements 
set forth by RWQCB and DHS regarding the use of recycled injection water would be 
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met. However, additional experiments should be conducted to increase the 
confidence in the results presented herein. There are several ways upon which the 
study may be improved. First, a one-day time step or smaller should be used in order 
minimize numerical dispersion errors. When using a smalt time step it may be useful 
to make use of the parallel version of FEM WATER in order to reduce run-times. 
Second. a detailed analysis of the effect of varying the head boundary conditions 
should be performed in order to provide the widest range of possible tracer travel 
times. Finally, the use of naturally occurring isotopic data to estimate groundwater 
travel times in the ABP may be useful in helping to validate future analyses. 
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